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Abstract. We study the resilience of complex networks against attacks in 
,,15 ' which nodes are targeted intelligently, but where disabling a node has a cost 

O 1' to the attacker which depends on its degree. Attackers have to meet these costs 

with limited resources, which constrains their actions. A network's integrity is 
^ , quantified in terms of the efficacy of the process that it supports. We calculate how 

w . the optimal attack strategy and the most attack-resistant network degree statistics 

, ' depend on the node removal cost function and the attack resources. The resilience 

(^ ' of networks against intelligent attacks is found to depend strongly on the node 

O ' removal cost function faced by the attacker. In particular, if node removal costs 

* J-»^ , increase sufficiently fast with the node degree, power law networks are found to 

^1,^. be more resilient than Poissonian ones, even against optimised intelligent attacks. 
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psj ■ 1- Introduction 

. -^ , In recent years there have been several studies into the resihence of complex networks 

(^ • against random failures and targeted attacks, in which a fraction of the nodes or of the 

f^ I bonds are removed. It was found that scale-free networks (with degree distributions 

that decay slowly via power laws, as in preferential attachment models) are more 

robust against random node removal than Poissonnian (or Erdos-Renyi) graphs, which 

may explain why many real-world complex systems involve networks with power-law 

. , distributed degrees. However, scale-free networks were found to be very vulnerable 

j^ ■ to intelligent attackers that target high-degree nodes [J [5]. Against edge removal, 

Poissonian and power law networks turned out to produce similar responses [3] . There 

are two reasons why we aim to study network resilience further. First, while the 

motivation behind such studies is that the networks provide the infrastructure for 

some process (with interacting 'agents' or processors occupying the nodes), and that 

process disruption is the true goal of an attacker, most authors measure the impact of 

attacks indirectly, via topological properties that serve as proxies for the integrity of 

the process (e.g. the overall connectivity and path- length statistics [4], or percolation 

characteristics [H El [5l [6j [7] ) . Here we seek to quantify the damage inflicted by attacks 

directly in terms of the process which the network is meant to support, similar to [8]. 
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This requires solving stochastic processes on complex networks with arbitrary degree 
distributions, which is what statistical mechanics enable us to do. Our second and 
most important reason is that network resilience has so far been studied strictly in 
the context of random or intelligent removal of a fixed fraction of sites or bonds. 
This seems unrealistic. In most real-world scenarios (attacks on computer networks, 
viruses attacking cellular networks, etc) attacking a highly connected node demands 
more effort on behalf of the attacker than removing a weakly connected one. Similarly, 
any sensible defender of a network would devote more resources to the protection of 
'hubs' than to the protection of 'outpost' nodes. The study of network resilience 
against attack or dilution calls for more appropriate and realistic definitions, that 
include the inevitable resource constraints faced by attackers and defenders alike. 

Turning to a formulation where attackers have finite resources, to be deployed 
intelligently when the cost of removing a network node depends on the degree of that 
node, changes the game drastically. It introduces a trade-off between the merit in 
terms of inflicted damage of targeting high-degree nodes versus the disadvantage of 
associated cost (attacking many 'hubs' may be inaffordable). One would like to know 
the maximum amount of damage that can be inflicted (by e.g. a virus to a biological 
network), given the limited resources available to the attacker (e.g. food, lifetime) and 
given the network's degree-dependent node removal costs. Similarly one would like to 
identify the most resilient network degree statistics to withstand an optimal attack. 
The answers to these questions may aid our understanding of structural properties 
of biological (e.g. proteomic) signalling networks, where competition and natural 
selection act as driving forces towards attack resistance, but also to aid the design of 
attack-resistant synthetic real-world (e.g. communication) networks. 

Here we develop a framework for the study of network resilience that includes 
limited attack resources, degree-dependent node removal costs, and resilience measures 
based on process integrity. We consider two types of processes where structurally 
different interacting variables are placed on the nodes of networks with arbitrary 
degree distributions: interacting Ising spins (where global order is ferromagnetic or of 
the spin-glass type) , and coupled Kuramoto oscillators (where global order is measured 
by synchronization). Both are solvable using finite connectivity replica theory, which 
enables us to quantify their integrity by the critical temperature of the ordered state. 
An attacker with finite resources seeks to destabilize these processes by removing or 
disrupting selected network nodes using his knowledge of the network's degrees. The 
attacker is also allowed to disable nodes partially (with a proportional reduction in 
attack costs). We identify the most damaging attack strategy, given a network's degree 
distribution and given the degree dependence of the node removal costs and the attack 
resources available. We then determine the optimal network topology from the point of 
view of the defender, i.e. that degree distribution for which the integrity of the process 
is preserved best when attacked by a foe who employs the most damaging attack 
strategy. The optimal attack strategy and the optimally attack-resistant network 
topology are found to be universal across the types of microscopic variables and types 
of global order considered. As expected, the resilience of network processes against 
intelligent attacks depends strongly on the node removal cost function faced by the 
attacker. Moreover, in sharp contrast to the traditional set-up where attackers are 
allowed to remove a fixed fraction of the nodes (and hence can simply target the 
'hubs'), we find that if node removal costs increase sufficiently fast with the node 
degree, and if attackers have finite resources to meet these costs, power law networks 
are more resilient than Poissonian ones, even against optimised intelligent attacks. 
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2. Definitions 

2.1. Processes, supporting networks, and constrained attack variables 

We study two systems in which interacting stochastic variables are placed on the N 
nodes of a complex network. The network is defined via variables Cij G {0, 1}, with 
Cij = 1 if and only if the nodes i and j are connected. We define Cy = Cji and ca = 
for all («, j), and abbreviate c = {cy}. The first system (A) consist of N Ising spins 
ai E {—1,1}, in thermal equilibrium, characterized by the following Hamiltonian 

A : H{cr) = - ^ Cy J^J^^^Jcr,aJ (1) 

i<j 

with <T = ((Ti, . . . ,(T]sr). The second system (B) consist of N Kuramoto oscillators, 
with phases 9i G [— tt, tt], again in equilibrium but now with the Hamiltonian 

B : H{e) ^-J2 c^3 'hM, cos(0, - Bj ) (2) 

i<j 

with = {di, . . . ,9n). The bonds Jy G IR are drawn randomly and independently 
from a distribution P{J). The variables ^i G S C [0, 1] in (|ll2p represent the impact 
of attacks, with ^^ = if node i is removed completely and ^i = 1 if it is left alone. We 
demand that 1 G S, so leaving a node intact is always an option, and for simplicity 
we take S to be discrete and finite. We define the node degrees fci(c) = ^ Cy , so the 
degree distribution and the average connectivity of c are p{k\c) = N~^ ^^ Sk,ki{c) a-nd 
(fc) = J2k>o f^PiMc)^ respectively. 

We assume that the cost to the attacker of setting ^,; = ^ at a node of degree 
ki — k is ip{^, k) > 0, where "0(1, k) = and dip{^, k)/d^ < for all ■^ G S. If attackers 
have limited resources they can only disrupt a subset of the nodes, since the {^i} will 
now be subject to a constraint of the form J^i i'i^ii h) < C. A natural choice for ip is 

^{^,k)^n{l-OHk) (3) 

where (l){k) is a non-decreasing function, with 0(0) = and 0(fc>O) > 0. The attack 
cost for a node increases with the number of links to/from it; disconnected nodes 
can be attacked for free. The normalization factor k is chosen such that the resource 
constraint takes the simple form iV"^ X^i V'(Ci: ^i) ^ 1- The attacker is assumed to 
act intelligently, using knowledge of the network's degrees, so the degrees {ki} and 
the attack variables {^i} will generally be correlated. Finally we draw the network c 
randomly from a maximum-entropy ensemble defined by a probability distribution in 
which the degrees are constrained to take prescribed values k = (fci, . . . , k^): 

Prob(c) = Z-i[k]n^fc,,fe,(c)> Z[k]=Y,l[h,Mic) (4) 

i c i 

We abbreviate p{£,,k) = N~^J2i^ii,i^ki,k, and define q{^\k) via p{£,,k) = 
q{£,\k)p{k). The resource constraint on the attack variables then translates into 
^cj.^(^, fc)(7(^|A;)p(fc) < 1. The attacker is assumed to know the degree sequence 
k of the network to be attacked, and can adapt accordingly the conditional likelihood 
g(^|fc) to maximize the impact of his actions; q{^\k) constitutes his attack strategy. 
The realistic regime is that where ■)/'(?, k) obeys J2k ^('^j ^)p(^) > 1' ^o that the trivial 
'destroy-air attack strategy q{£,\k) = ^^.o Vfc is not feasible (i.e. too costly). 
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2.2. Quantifying process integrity and optimal attack and defence strategies 

With each process (|ll2p running on the network and each associated ordered phase 
(ferromagnetic, spin-glass, or synchronized) corresponds a critical temperature Tc, 
which will for large N depend on the network and attack characteristics k and ^ only 
via q{^\k) and p{k). The larger Tc, the more robust is the ordered phase against 
local noise, so we can quantify the integrity of the process by the value of Tc[p, q]. The 
attacker wants to destroy the ordered phase of the process, whereas the defender seeks 
to protect it. This allows us to give precise definitions for the optimal attack strategy 
and the optimally resistant degree distribution in terms of process integrity. The 
optimal attack strategy q*[p] is the conditional distribution q{£,\k) for which Tc[p,q] is 
minimal, given the degree distribution p and given the resource constraint: 

q*[p] - argmin^^^^^^ ^(e,fc),(«|fe)p(fc)<l}^c[p, q] (5) 

The optimal (most resistant) degree distribution p* to be chosen by the defender, 
given the average connectivity c (finite network resources) and attack cost function ^ 
is then that p{k) which subsequently maximizes this g-niinimized critical temperature: 

p* = argmax^p^^^ kp{k)=c}Tc[p, q*[p]] 

= argmax^^^^^ ,p(,)^^jargmin^^_^^^ ^(e.fe).(e|fe)p(fe)<i}^-[P' ?] 

(6) 

The end result is a situation where the defender, by choosing an appropriate degree 
distribution, maintains the highest achievable critical temperature Tc[p* , q*[p*]], given 
he is subjected to the most damaging attack. However, within this scenario one could 
in fact ask many more interesting questions, such as what would be the effect of 
misinformation, a situation where a defender optimizes the network on the basis of an 
anticipated attack q*[p] (so he chooses degree distribution p*) but is then faced with 
an attack with strategy q' ^ q*[p], so that the actual critical temperature is Tc[p*, q']. 
We see that the problem of identifying the optimal attack and defense strategies 
(|5l6p splits automatically into two distinct parts. The first part is calculating the 
critical temperature(s) Tc[p,q] of the relevant phases. This is done by evaluating for 
the systems ()1I2|) the asymptotic disorder-averaged free energy per spin /, from which 
one extracts the phase diagrams for systems on typical graphs from (j4|: 

7a^ - lim ITTT^ogY e-Pfii'^) (7) 



Jb= - lim ^log / d0 e-/5^(e) (8) 

in which (3 — T^^ (where T denotes the temperature), and where "^^^ denotes averaging 
over the disorder in the problem, viz. the randomly drawn graphs with statistics (H]) 
and the random bonds {Jij}. The calculation of ([7]) and ([5]) is done with the finite 
connectivity replica method, based on the identity logZ — lim„_j.o n^^ logZ", and 
details are relegated to [Appendix A| and [Appendix B| in order not to disrupt the flow 
of the paper. The second part of the problem, to be tackled once the formulae for 
Tc [p, q] have been derived (which, expectedly and fortunately, turn out to be simple 
and very similar across models and ordered phases), is to carry out the constrained 
optimizations in (J5l6p . by a combination of analytical and numerical techniques. 
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3. The process integrity measure 

We show in the appendices of this paper that the critical temperatures Tc[p,q] for 
the emergence of global (F or SG) order, given we choose the bond distribution 
P{J) = l{l + r])S{J-Jo) + i(l-77)(5(J+ Jo) (with Jo > 0), follow for both Ising 
spins and coupled oscillators from formulae of the following form: 

F: aW,(/3) = 1, A(i): eigenvalues of mW(/3) = r;i^(/3JoCC')7(e') (9) 

SG: X[ll{f3) = l, A(2): eigenvalues of Afg,^(/3) = i^^^^j^^^,)^^^,) ^^0) 

in which /3 = 1/T and 

k 

Here K{z) = tanh(z) for interacting Ising spins and K{z) = Ii{z) / Iq{z) for coupled 
oscillators. In both cases K{—z) = ~'K{z), ^K{z) > 0, and liuiz^oa K{z) = 1. 
There is no F phase if r; < 0, so we take rj > from now on. The structure of the 
above formulae is in agreement with results from percolation theory and spreading 
phenomena, which show that the threshold characterizing the percolation transition 
or an epidemic outbreak in a network depends on the ratio {k'^) /{k) of the first two 
moments of its degree distribution [TJ O [6l [3 [2j [9l [TOl [11] . The approach followed 
here is closer to the envisaged picture of interacting agents or processors on network 
nodes, and has the benefit of applying to the whole interval S = [0, 1], as opposed to 
S — {0, 1} which can be accessed by percolation theory. 

3.1. Tests and bounds for critical temperatures 

Before any attack one has ^ G {!}, so 7(^) = 7(1) — (fc^)/(fc) — 1 and the above 
formulae would have reproduced the known results for the unperturbed system, viz. 

F: rjK{l3Jo)[{e}/{k}-l] = l (12) 

SG: K^(3Jo)[{e)/{k)-l] = l (13) 

Another simple test is to consider S — {0, 1}. Here each node is either unaffected 
or removed completely, leaving a new network identical to an unperturbed network 
as described by (|12|13|) . but with reduced size N' = J2i^i' ^^'^ '^ith degrees 
K — X^j Cyfj. We would find, in the case of random attacks q{£,\k) = C'5c,o + (l^C)'5c,i- 

(fc)' = li- TT^E^'--^. - (i-O(fc) (14) 



N^oo {1~C)N 



ij 



{ky - ltn^-J—-J2^^c,,ca^,^i = (1-C)'(fc2) +C(l-C)(fc) (15) 

giving the following transparent formulae for the post-attack transition points: 

F: ,7if(/3Jo)(l-C)[(A;')/(fc)-l] = l (16) 

SG: K^(3Jo)il-0[{k^)/{k)-l] = l (17) 

If, alternatively, we apply to this scenario the result (|9ll0p . we find 7(^) = 
q{^)[{k^) / (k) —1] and K{13Jq^^') — K{f3JQ)6^^i5^i^i, and the relevant matrices reduce to 
M«(/3) =. rjKif3Jo){l-C)[{k^)/{k)-m,iS^,', andMg)(/3) = KHf3Jo)il-Cm')/{k)- 
l]'^€,i'^C',i- 0^6 solves the eigenvalue problems trivially, and indeed recovers (fT6lfT7|) . 
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A final trivial test is to consider q{^\k) — 5^^^^, where ^o G (0; 1)? ^n attack equivalent 
to replacing Jq -^ QJo- Upon substituting this choice into (|9ll0p one confirms, via 
7(C) = <^{,{o[(^^)/(^)~l]' tti3,t our general theory indeed reduces to (|12ll3p with the 
correctly reduced coupling strength. 

Solving the eigenvalue problems (|9I10I) analytically is not always possible, but 
eigenvalue bounds are obtained easily. Our matrices are of the form M^^> = 
L(^^')7(^'), where L{u) — riK{/3Jou) for the F transition (so L{u) is anti-symmetric) 
and L{u) = K^{(3Jou) for the SG transition (so L{u) is symmetric), and where 7(^) > 
for all £,. We symmetrize the eigenvalue problem Xx{£,) = ^,, M^^>x{^') by defining 

2/(0 = ^iOVW), giving XyiO = Ei'[VliOL{^^')VW)m')- This implies that 
Amax = max =^ — r— T (18) 

y L^r(t) 

which can be simplified to 



E^^'y{OVW)me\)VW)yie) 



Amax = max — ^S ^ — (19) 

y E^y (0 

Variational arguments can now be applied in order to get lower bounds. In particular, 
upon substituting j/(^) = 6^ f and varying ^ one derives the statement 

^ {lione)} (20) 

To find upper bounds, we use the fact that the maximum in (|19p will have y(^) > for 
all ^. We then use the inequalities L{\u\) < ar]PJo\u\ (for F) and L{\u\) < (a/3Jo)^|up 
(for SG), where a = 1 for Ising spins and a = | for coupled oscillators, to get 

KL < (^VpJo max <^ ^ = — \ (21) 

.(,) ^. ,,,2 r E,ey(07^(oienfP7^(e)2/(f) i _, 

ALL < (apJo) max <^ ^ = — j-- \ 22 

The last two maxima are calculated easily, leading us to 

F: T,[p,q] < 77aJo^e'7(0 (23) 

? 

SG: T,[p,q] < aJo(^eS(C))^ (24) 

3.2. Explicit simple form for a process integrity measure T[p, q] 

The inequalities p3l24p become equalities for large c, where the critical temperatures 
diverge and hence /3 — >• in (fT2| [T3\\ ; the right-hand sides of (|23|24p then become the 
true integrity measures of the process. Moreover, for certain natural choices of the set 
S the latter statement is in fact true for any connectivity c. For instance, if S C {0, 1} 
(all nodes are either fully disabled or left alone) one may use if (/SJo^O — S,S,'K{(3Jq) 
to diagonalize the matrices in (I9ll0p and find 
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which reveals that the critical temperatures are monotonically increasing functions of 
the sums ^. S.'^liO ^^^ F-type order and J^e C^7(C) ^oi SG-type order (for S — {0, 1} 
the two sums are in fact identical). In view of these properties, and in view of the 
minor differences between the F and SG cases, in the remainder of this study we adopt 
the quantity J^f C^7(C) ^^ our integrity measure, giving 



nP,<l]^TjT}_^^'QimPik)k{k-l) (27) 

We define the set of relevant degrees k a,s S = {k > 1| p{k) > 0}. The optimal attack 
strategy is then the choice q*[p\ which solves the following optimization problem: 

minimize: r[p, (/] (28) 






subject to : q{^\k) > V(^, k), }^ q{^\k) = iykeS (29) 

E J2(^-0^ik)qimP{k) < ^-' (30) 

CeHfeeS 

To avoid trivial pathologies we assume that 3k > 2 with p{k) > (if untrue we would 
not have an ordered state in the first place, as it would have given Tc[p, q] — 0), and 
that 'YlikPik-)k{k — 1) < oo (if untrue there would not be a finite critical temperature 
before the attack). Clearly q*{S,\k) = 5^i_ for k ^ S; any other choice would sacrifice 
attack resources without benefit. The best defense against optimal attacks is the choice 
for the degree distribution p{k) such that the above minimum over q is maximized. 

3.3. Bounds on the process integrity measure 

To judge the quality of attack strategies it will prove useful to have bounds on 
the value r[p, g*[p]] corresponding to the optimal attack q*[p]. An upper bound is 
easily obtained by inspecting the result of non-intelligent random attacks of the type 
q{^\k) = (1 - Q)S^^o + QSi,i, with < Q < 1: 

r[p,q]^[{k^)/{k)-l]Q (31) 

> 1 - l/K{cj){k)) (32) 

The sharpest bound of this form follows when seeking equality in the last line, givingj: 

r[p,q*[p]] < r[p,gra„dom i] = (l - ^(^)<^(^-l)>/<^> (33) 

If S = {0, 1} then p3p is the best possible upper bound based on random attacks. 
If S = [0,1] we can improve upon ([55]) by investigating random attacks of the form 
9(^1^) ^ '^K ^ C]- The optimal choice turns out to be ^ = 1 — 1/k((/)), giving 

r[p,q*[p]] < r[p,<7;,„dom 11] = (1 - ^^(^)V^(fc-i)>/(fc> (34) 

To find lower bounds for r[p, q*[p]] we first define modified probabilities 7r(/c) G [0, 1]: 

^, , _ EggE(i~e)g(e|fc)0(fc)p(fc) ,„_, 

"^ ^" {E,,di-omk)m) ^ ^ 

X Note that K{<f>{k)) > 1 due to our earlier ruling out of the trivial attack strategy q{£,\k) = <5^,o. 



> 
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with associated averages written as (. . .)^. Note that the denominator of 
bounded from above by k~^, via the resource constraint. We can now write 

r[p,q] =. ^M^ ^ ±.Y.p(^k)k{k-i)^{i^eh{m 

(fc(fc-l)) 1 / fc(fc~l)Egg5(l-^')g(C|fc) \ 

(fc) K{k)\ 0(fc) Eees(l-09(ei^) /- 

^ (fc(fc-l)) CS / fc(fc-l) \ .^.^ 

(fc) k(/c) \ 0(fc) /.r ^ ' 

in which the factor Ch > depends only on the choice made for the value set S: 

CH=max|i-^|, with (/(O).. =E^(^)/(^) and ^^0 = 1 

(37) 

One easily proves using S C [0, 1] that Ch £ [I, 2], that C{o.i} — l^ ^nd that C[o ij = 2. 
We conclude, in combination with p3l34p , that 

^^rnn- 1 1 / fc(fc-i) \ ^ (fc)r[p,g-[p]] ^ ^ 

" ■^ ' -^' K(fc(fc-1)>\ (/)(fc) At- (fc(fc-l)) - K(<^(fc)) ^ '^ 

^^foii- 1 2 / fc(fc-i) \ < (fc)rb,g*b]] < A ^f ..q. 

" ^' ^' K(fc(fc-1))\ (/)(fc) L- {k{k-l)) - V n{(t){k))) ^ ' 

The lower bounds are satisfied with equality if the attack resources are exhausted and 
if ((1 — ^))q — ((1 — 0)5 for each q(i,\k) with k G S; the last condition is always met 
if S = {0, 1}. However, the lower bounds still depend on the attack strategy via the 
measure tt. From p8|39p and the general property r[p, q] > 0, which follows from the 
definition of T[p, q], we finally obtain the strategy- independent bounds 

H.,oa|: ,„..{„, :-^} , <|iiM , :--^ ,0) 

, , r 2i?/K ^ (k)r\p,q*\p]] ( 1 \2 , , 

5=0,1 : maxiO, 1- ,, ,, ' ,, \ < ^ ' ,, ' ^ ,''" < (l- ,,,,,J 41) 

^ ' ^ I ' (fc(fc-l)>/ - (fc(fc-l)> ~ V n{(t>{k))J ^ ' 

with 

i? = maxfce5{fc(fc-l)Mfc)} (42) 

The latter bounds reveal immediately two distinct situations where it is not possible 
for any intelligent attack to improve on the damage done by random attacks: the case 
4>{k) = k{k — 1) for all fc S S* (here the benefit of degree knowledge exactly balances 
the cost to the attacker of using it), and the case of regular random graphs, viz. 
p{k) — (5j, /fc\, where there is no degree knowledge to be exploited in the first place. 

4. Optimal attack and optimal defense for S — {0, 1} 

The attacker's objective is to minimize r[p,q]. We have seen that for 5 — {0,1}, 
where nodes are either fully disabled or left alone and C'lo,!} = 1, the lower bound in 
(|40p could in principle be realized. This will serve as an efficient guide in finding q*[p]. 
Attack strategies for S = {0, 1} are of the form q{^\k) = q{0\k)S^fi + [1 ^ 'Z(0|fc)]'5^,i, 
so we need to determine q{0\k) for all k £ S. 
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4.1. Construction of the optimal attack strategy 

We first define the attacker's 'target' degree set A'^ S, witli R as defined in P^ : 

A^{keS\k{k~l)/(l){k)^R} (43) 

Tlie inequality (fc(fc— l)/0(fc))7r < R used in tfie final step of our derivation of (|40|) is 
satisfied with equality only if n{k) — for all k ^ A. According to (1551) this requires 
q{^\k) = S^^i for all k ^ A. The only remaining requirement for satisfying the lower 
bound in pO)) is that we satisfy the resource constraint with equality. Hence the set 
of optimal attack strategies is defined strictly by the following demands: 

\JkiA: q{0\k) = (44) 

\/keA: g(0|fc)e [0,1], ^(j(O|fc)0(fc)p(fc) = l/K (45) 

keA 

It is straightforward to verify directly, using (pik) = k{k — l)/R for all k ^ A, 
that strategies satisfying these conditions indeed give the lowest possible value for 
T[p, q] according to our bounds, and satisfy the resource constraint with equality. By 
construction, the set A cannot be empty. 

At this stage in our argument we must distinguish between two distinct cases. In 
the first case the attacker need not look beyond nodes in the target set A (031), since 
removing those will already exhaust or exceed his resources; he will simply remove as 
many of those as can be afforded. In the second case the removal of all nodes in A 
does not exhaust the attack resources, and new target sets need to be identified: 

• The target set A is exhausting, X^fce^t 4'{k)p{k) > 1/k: 

Here it is immediately clear that optimal attacks will indeed exist, i.e. the 
conditions (|44|45|) can be met. Only nodes from A will be removed. If there 
is at least one k* E A with (t>{k*)p{k*) > 1/k, the attacker can simply execute 

k* = argmaxj.g^{0(fc)p(fc)} (46) 

q{0\k*) = 1/Kp{k*)(t){k*), yk^k* : qiO\k) = (47) 

If instead (j){k)p{k) < 1/k for all k & A there is no target degree in A which 
would on its own exhaust the attacker's resources. The attacker will first remove 
all nodes with degree fcjf = a.rgniax^^^{(j){k)p{k)} by setting q{0\ki) = 1. He 
will next direct attention to the reduced set A/{ki} and remove nodes with 
degree ^2 = ^^E^^^k£A/{k*}{4>{k)pik)}, etc, until the resources are exhausted. 
At the end of this iterative process the attacker will have removed a sequence 
of degrees {k*, • • • , k^} C A (where nodes with degree fc^ will generally be only 
partially removed, as allowed by remaining resources). In words: the attacker 
first determines the target set A of those degrees with p{k) > for which the 
ratio k{k — l)/(f>{k) is maximal. He then ranks the degrees in A according to 
the value of 4>(k)p(k) and proceeds to remove degrees iteratively according to 
this ranking until his resources are exhausted. This strategy will always lead to 
{k{k — \) / (f){k)) T,: = i?, and satisfy the lower bound in (PU)) with equality. 

• The target set A is non-exhausting, X^fee^ 4'{k)p{k) < 1/k: 

Here the attacker can afford to remove completely all degrees in the set A, but 
setting q{0\k) = 1 for all k G A does not exhaust his resources. He should 
subsequently direct attention to those nodes in the reduced set S/A for which the 
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ratio k(k — l)/0(fc) is maximal, and so on. The result is again an iteration, at the 
end of which the attacker will have removed a set of degrees {fc*, • • • , fc^} D A 
(where nodes with degree fc^ will generally be only partially removed). In this 
case (fc(fc — \)/4>{k)) < R and the lower bound in (PH)) is no longer satisfied 
with inequality; however, this does not imply that the strategy is non-optimal, 
since it might be that the bound is no longer tight. Here it is therefore difficult 
to prove rigorously that the identified strategy always constitutes the optimal 
attack, but it is the logical continuation of the optimal attack identified earlier 
and its optimality is consistent with numerical experiments (to be shown later). 

We can combine both cases above in a transparent iterative attack protocol. We define 
at each step t. the target set Ai, the set Si of nodes that have not yet been targeted, 
and the resource remainder A^ — k~^ — ^^ q{0\k)ip{k)p{k) . The process is initialized 
according to S'o = S* and Aq — k~^, and starts with q{0\k) = for all k. It is iterated 
until Af = 0, according to 

step 1: calculate new ratio Ri = maxfcg5^_j{fc(fc— l)/(/)(fc)} 

step 2: identify target set Ae — {k G Si-i\ k{k—l)/(j){k) — Rg} 

step 3: choose (any) fc^ e At for which (j){k*g)p{kl) = maxfc£^^{0(fc)p(fc)} 

step 4-' check whether attack resources can be exhausted: 

<t){k*g)p{k*g) > Af_i : yes, 

remove as many degree fc^ nodes as possible 

set q{0\kt) = Ai.i /(j){kt)p{k*g) 
A^ = 0, attack terminates 

4,{kl)p{k}) < A£_i : no, 

remove all degree k^ nodes 
set q{0\k*i) = 1 
step 5: define Si = Si-i/k*^ and A^ = Af_i — 4>{k*g)p{k*i)q{Q\k*g) 

Always the end result is a sequence {k\^ . . . , /c^.j^} of target degrees that are fully 
removed, possibly supplemented by a further degree fcj, of which a fraction will be 
removed (to exhaust fully the attack resources). 

^.2. Properties of the optimal attack strategy 

We next evaluate the impact of the above attack strategy q* \p\ on our process integrity 
measure. We define the set A* = {k\, . . . , fc^_i} of fully removed degrees, and write 
k*j^ simply as k* . The post-attack value of the process integrity measure will be 

r[p,g^H] = ^^^^y^ " (^ ^ p{k)k{k-i) ±.q[Q\k*)k\k*-i)p{k*) 
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The attack q*[p] exhausts all resources, so Ai_ 



T[p,q*[p]] = 



(fc(fc-l)) 
(fc) 



(k) ^ 



Pik)cb{k) 



Efce^* '/'(fc)p(fc)- Hence 
fc(fc-l) k*{k*-l 



<t>{k) 



(/){k*) 



k*ik*-l) 



K{k)<i){k*) ^'^^^ 

Since by definition k{k-l)/(j){k) > k*{k* -l)/4>{k*) for ah k e A*, both the second 
and the third term of P^ are strictly non-positive. 

The result (HHl) can be compared to that of the most damaging random attack 
(where no degree information is used), viz. to (p3)) . The benefit AT = r[p, g*[p]] — 
'^[Pj ^random] ^^ ^^"^ attackcr of using optimal intelligent attacks as opposed to optimal 
random attacks then takes the form 



Ar 



-- V p{k)c^{k) 



1 

(fc) 

1 

(fc) 



E 

keA 



E 

keA* 



— V p{k)cp{k) 



\k{k~l) k*{k*-l) 
(j){k) 0(fc*) 

k(k-l) k*{k*-l) 



1 



k*{k*~l) (fc(fc-l)) 



K{k) I (j){k*) 



Wk)) 



0(fc) 



1 



K{k) (0(fc)) 



0(fc*) 

J2p{k)^{k) 



1 



E 

keA* 



p{k)cj,{k) 



fc(fc-l) fc*(fc-^-l) 
</.(fc) " 0(fc*) 

~ p(fc)0(fc) 



1- 



fc*(fc*-l) A:(A:-1) 
0(fc*) 

1 



4>{k) 



K{k) (0(fc)) 



fc^yt- 



«;(0(fc)) 

k{k-l) k*{k*-l) 



0(fc) 



0(fc*) 



(50) 



Since the set A* C S* is constructed specifically from those degrees for which 
k{k ~ l)/(j){k) is maximal, and since (0(fc)) > k~^, both terms of AT are strictly 
non-positive. One will thus generally have AF < 0. Again we also recognize the two 
special cases where there will be no gain in intelligent attacks, namely 0(fc) = k{k — 1) 
(with any degree distribution), and p{k) — Sk,{k) (with any cost function (f>{k)). 

In terms of the dependence of our results on the cost function (j){k) it is clear that 
everything evolves around the dependence on k of the ratio (j){k)/k{k — 1). This ratio 
represents for each k the balance between the cost of removing degree-fc nodes versus 
the benefits in terms of damage achieved. If for simplicity we choose (f>{k) = k'^{k~l), 
then for C < 1 the intelligent attack will be to take out first the nodes with the largest 
degrees k G S that can be removed without violating the resource constraint (i.e. the 
'greedy' attack strategy is optimal), whereas for ^ > 1 the intelligent attack will target 
first the nodes with the smallest degrees k G S that can be removed without violating 
the constraint (here attacking hubs is too expensive to be efficient). Furthermore, it is 
not at all a priori clear what would be the most resistant degree distribution against 
such attacks in the presence of resource constraints. Naively one could perhaps have 
expected that for small C (where the attacker will target hubs) the best strategy for 
the defender could be to choose a narrowly distributed degree distribution, so there 
are no hubs to be exploited. Interestingly, we will see below that that is not the case, 
and the optimal degree distribution can be more subtle. 
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Species 



N 



(fc> fcr 



Method Reference 



C. Elegans 


3512 


3.72 


524 


Y2H 


il2j 


C. Jejuni 


1324 


17.52 


207 


Y2H 


m 


E. Coli 


2457 


7.05 


641 


PMS 


m 


H. Sapiens 


9306 


7.53 


247 


HPRD 


mi 


S. Cerevisiae 


3241 


2.69 


279 


Y2H 


TB" 



Table 1. Size A^, average connectivity {k), maximum degree fcmax, experimental 
detection method or source, and reference for the biological (protei interaction) 
network data sets used in our numerical experiments. The detection methods or 
sources are abbreviated as follows: Y2H, Yeast two-hybrid; PMS, Purification- 
Mass spectrometry; HPRD, the human protein reference database. 



5. Numerical results 



5.1. General methods 

In this section we ihustrate, apply and extend via numerical experimentation the 
results derived above. We determine by numerical maximization the most resistant 
degree distribution against optimal intelligent attacks, and we compare for typical 
biological networks the effects of optimal intelligent attacks in terms of process 
integrity against optimized random attacks and against the bounds established earlier. 
The biological networks used are experimentally determined protein interaction 
networks (PINs) of different species, namely C. Elegans, C. Jejuni, E. Coli, H. Sapiens 
and S. Cerevisiae; see Table [T] for characteristics and references. For each biological 
network we also generate several synthetic alternatives with the same size N and 
average connectivity (fc) as the biological one, but with different degree distributions: 
Poissonian, the optimally resistant degree distribution, or a distribution generated 
via preferential attachment with a fat tail similar to the biological network. In 
all cases we choose node attack cost functions of the form (j){k) = k'^{k — 1), with 
C = 0, 1, 2, and we set the attack resource limit to k~^ = {ip{k))'p/q, where g > 1 is a 
control parameter and the average (. . .)-p is calculated over the Poissonian distribution 
V{k) — (fc) e^^'^ /k\ For each PIN and each synthetically generated counterpart the 
average cost function {4>{k)) is found to be always larger than or equal to the one 
calculated over the equivalent Poissonian distribution (if C = 0, where (/>(fc) = fc — 1, 
equality of course holds trivially for all distributions since they share by construction 
the value of (fc)). This ensures that for q > 1 the attacker's resources will in all our 
experiments be in the relevant regime (0(fc)) > k~^. The optimally attack resistant 
networks are found via a stochastic graph dynamics, starting from a biological protein 
interaction network, in which at each step a bond is selected at random and is moved 
to another location if this move increases the post-attack integrity measure r[p, (Z*[p]]. 
Bond relocations are the minimal moves that preserve the average degree of the 
network. After each move, the optimal attack strategy q*[p\ defined in the previous 
section is applied to the new network. In order to prevent the graph dynamics from 
getting stuck in suboptimal configurations, we allow initially for groups of bonds to 
be moved, and as the algorithm proceeds the size of these groups is reduced, in the 
spirit of |T7]. 
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Figure 1. Left: log-log plot of the degree distribution of the H. Sapiens PIN, 
synthetically generated networks with the same size A'^ and average connectivity 
(k) as H. Sapiens PIN but with Poissonian and preferential attachment degree 
distributions. Right: log-log plot of the degree distribution of the network that 
has the same size and average connectivity as the H. Sapiens PIN, but that has 
been constructed to be optimally resistant against optimal intelligent attacks, 
given the node removal cost function <f>{k) = k'^{k — 1), with <^ = 0, 1, 2, and given 
available attack resources, characterized by re^^ = -{<l>{k))-p, with q = 3, 5 (see 
legend). We observe that upon decreasing (^, where the optimal attack strategy 
starts targeting the high-degree nodes, the optimally resistant degree distribution 
changes takes a binary form, describing a module of high-degree nodes in a sea of 
unconnected nodes. 



5.2. Degree statistics before and after network optimization 

Figure [T] shows the resuhs of applying the above procedures to the H. Sapiens PIN, 
for C = 0, 1, 2 and q = 3,5. For C = 1, 2 (where it is not advantageous to the attacker 
to target high-degree nodes) the optimally resistant degree distribution p* [q] is seen to 
exhibit a smooth dependence on the degree k. For C = 0, where the degree dependence 
of node removal costs is modest and the optimal attack strategy is to target high degree 
nodes, one could expect the optimal network to become regular, in order to disallow 
attackers to benefit from degree information. Instead, we observe an entirely different 
solution. Here, the optimal defender produces as many hubs as possible, so that the 
attacker is unable to remove all of these. The result is a distribution of the form 



p*ik) 






K 



5k. 



K, 



K> (k) 



(51) 



with K = 192 for a network with size and average connectivity identical to the H. 
Sapiens PIN. In a situation where attackers can and will target nodes with maximal 
degree first, it appears that the optimal defender chooses a network with a "modular" 
configuration, with a core of nodes highly connected to each other, in a sea of 
disconnected nodes. The attacker is prevented by resource limitations from removing 
more than a (tiny) fraction of the core. The strategy of the optimal defender is to 
sacrifice a few highly connected nodes to save many. 

For distributions of the form (fSTj) the attack resources can be exhausted, and the 
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optimal attack is q*{Q\k) — Sk.K K / k(J){K) (k) . This results in 

r[p,q*[p]]= K -l-K^-'^/K{k} (52) 

Insertion into (j40|) shows that both bounds are now satisfied with equality. For 
C > the defender would wish to choose K as large as possible, but for a finite 
network there is a limit. There are just N{k)/K connected nodes; if each of these 
is to have K neighbours we must demand N{k)/K — 1 > K, i.e. K < Kc = 
\jN{k) + 0{N^). The same result follows from general entropic arguments [TSl 119) . 
For large N the number of graphs with degree distribution p{k) equals exp[7VS'] where 
5=i(fc)[log(A^/(fc)) + l]-^^p(fc)log[p(fc)/7r(fc)],with^(fc) = e-<fc>(fc)Vfc!- For dSH 
one obtains 

S^]^{k)[l-\og[-^)]+0{K-HogK) (53) 



Again we obtain the cut-off point K < K^ ~ ^jN{k) for graphs with ([5T|) to exist. 
The value K = 192 found numerically, see Fig. [TJ is consistent with this bound (for 
H. Sapiens one has K^ = 264), but not identical to it. This is expected to reflect finite 
size corrections to our theory, and the fact that the theory requires all relevant k to 
be finite relative to N ^ whereas close to Kc one has k — Oi^/N). 

Note, however, that the distribution chosen by the optimal defender is not always 
exactly of the form (1511) . In some situations (depending on the amount of resources 
available to the attacker) the peak at fc = is not strictly 5-shaped, so that K is 
no longer subjected to the previously identified cut-off Kc, and one indeed observes 
the second peak to move to higher values of K (albeit with a reduced height). This 
results in a bimodal distribution with a (5-peak at some K > Kc and a broader peak 
at fc = 0, corresponding to a strongly disassortative network configuration (for the 
notion of assortativity see e.g. |5D]) where a small number of hubs are connected with 
an extremely large number of low degree nodes (reminescent of results derived in [21'). 
For this later distribution, as was the case with the bimodal distribution (|5ip with two 
strictly (5-shaped peaks, the attacker will again exhaust his resources upon removal of 
just a tiny fraction, q{0\K) — 1/k(J){K)p{K), of hubs. 

The actual distribution p* (fc) selected by the optimal defender when the attacker 
is bound, by resource limitations, to play the strategy g*(0|fc) = 6k, k ^ / i^4'{k)p{k) , is 
the one which maximizes the minimal integrity measure 



r(p,g*W) = — ( {k{k - 1)) - ^p(fc)fc(fc- i)g(0|fc) (54) 



^ ■\{k{k-i))-i 

\ k / 

{{k{k-l))-K-^K^~'^) (55) 



I ( fCi K" — 111 — 

(fc) 

1 

-(fc) 

achieved by the attacker. It is clear that the shape of the optimally resistant 
distribution will depend on the interplay between (fc^) and K, which is controlled 
by the resource limit k~^. For C — numerical studies show that for q sufficiently 
small (large amount of resources) p*(fc) assumes the shape (|5T|) . whereas for large q 
(small amount of resources) the width of the peak at fc = increases and the second 
peak moves to K ^ Kc- 

5.3. Values of process integrity measures before and after attacks 

In Figure[5]we plot the integrity measure T[p, q], for the different network distributions 
considered, before and after an optimal intelligent attack. We also show the values 
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Figure 2. Values of the process integrity measure r[p, 5] before attack (A), after 
optimal intelligent attacks (*), and after optimal random attacks (O), for different 
networks. The specific networks considered are experimentally determined 
PINs of different species (C. Elegans, C. Jejuni, E. Coli, H. Sapiens and S. 
Cerevisiae) and their synthetically generated counterparts with the same size and 
average connectivity, but different degree distributions (Poissoninan, preferential 
attachment and optimally attack resistant degree distribution following the attack 
q*lp])- The node attack cost function is </<(fc) = k''(k — 1) and the available 
attack resources are characterized by k~^ = {(f>{k))-p/q, with (^ = 0,1,2 and 
g = 3, 5 (see legends). The theoretical upper and lower bounds l|40p are shown 
as dotted and dashed lines, respectively. All results consistently reproduce the 
built-in order Fbcfore > Trandom > rintelligent (viz. A>0>*). Furthermore, the 
network realizations are consistently ranked, with the optimally resistant network 
(as expected) always outperforming the others, but with also the biological and 
preferential attachment network outperforming their Poissonnian counterparts. In 
fact, the degree of resistance of the optimally resistant network is quite remarkable. 
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that the integrity measure would take after a random attack (where sites are picked up 
at random and removed until resources are exhausted), with the same attack resource 
limit. The node attack cost function and resource limit chosen are 0(fc) = k''{k—l) 
and K~^ = {4){k))-p/q, respectively, with C = 0, 1,2 and q> 1. In addition we show 
the lower and upper bounds (I40p on the process integrity measure after an optimal 
intelligent attack q* \p] , as described by the protocol in the previous section, as dashed 
and dotted lines. As expected, the data points for random attacks always coincide 
with the dotted line of the theoretical upper bound (since such strategies formed the 
basis from which the upper bound was derived). Optimally attack resistant degree 
distributions are expected and indeed seen to be the ones for which the integrity 
measure of the network process after optimal intelligent attack is highest, but when 
quantified via T[p,q\ as in the figure one is struck by how well they perform, i.e. by 
the remarkably small reduction in the process integrity measure which they exhibit. 
For C, — 1^ the theoretical lower and upper bounds coincide with each other, and with 
the integrity measure values for random attacks and and optimal intelligent attacks. 
Here 0(fc) = k{k — 1), so costs and benefit for the attacker of degree knowledge 
balance each other out, and we have already showed that there is then no scope for 
the intelligent attacker to improve on the damage inflicted by random attacks. One can 
often understand the actual values obtained for T[p,q\. The relative reduction of the 
integrity measure before and after random attacks, for instance, can be calculated from 



before 



(|40]) . and for our resource limit k^^ = {(/){k))-p/q, this gives (rbcforc — rafter)/r 
{(j){k))-p/q{(j){k)). For C = this is always equal to 1/q; for ^ — 1,2 it is small 
for degree distributions with large second and third moments; much larger than for a 
Poissonian distribution. Similarly, for optimal intelligent attacks equation (I40p yield an 
upper bound on the relative attack-induced reduction of the process integrity measure: 

(rbefore-rafter)/rbeforo< [(/c? (/c - l))p/g(fc(fc - 1))] maXfeesIA^^"^}. For C = 1, 2, this 

change is again small for optimally resistant and biological networks, as a result of 
their large degree variance (except for C. Jejuni, which is distinct due to an unusually 
large average connectivity). 

Our results re-confirm that random and hub-targeted attacks have similar effects 
on Poissonian graphs, as often remarked in literature, due to the large homogeneity 
of the degrees. For regular graphs they would have produced identical results. 
However, one should be careful in concluding from this that processes running on 
Poissonian networks are hence the most resistant ones against hub removal. On 
the contrary, Fig. [2] shows that they are the most vulnerable ones, as their post- 
attack integrity measure is smallest. Interestingly, we find that processes running 
on networks produced by a preferential attachment mechanism are more resistant 
than those running on Poissonian networks, against both random attacks and optimal 
intelligent attacks. All this is due to the profound impact of resource constraints 
on the network resilience problem. Moreover, the degree distributions found in 
biological PINs generally exhibit, in turn, higher values for the post-attack process 
integrity measure than both Poissonian and preferential attachment networks, for 
random attacks and optimal intelligent attacks. A final feature emerging from Fig. 
[2] is that, while overall more robust compared to their preferential attachment and 
Poisonnian counterparts, biological networks seem significantly more resilient against 
random attacks than against hub-targeted attacks (see the top two panels with C = 0, 
where the optimal attack indeed targets hubs). 



Network resilience against attacks constrained by node removal cost 



17 



0.8- 

Ar 
r 0.6- 



0.4- 

o 

o 

0.2L- 






optimal 


A 


biological 


D 


preferential 


* 


poisson 



* * 



o O 



0.1 



0.2 
/ 



0.3 



0.8- 



0.6 



0.4- 



A^-^0 



00 

00 



0.1 



0.2 
/ 








optimal 


A 


biological 


D 


prcrcrenlial 


* 


poisson 




theory 



« **. 



0.3 



0.4 



Figure 3. Scatter plots of the relative variation AF/F of the integrity measure 
versus the fraction / of sites removed under optimal intelligent (left) and random 
(right) attacks, with 4>{k) = k — 1 and k~^ = (<j>{k))/3. Different markers 
correspond to different network families. The biological family is composed of 
the five experimentally determined PINs of Table [T] The other families are the 
synthetically generated counterparts of the biological PINs, with the same size and 
average connectivity, but different degree distributions (Poissonian, preferential 
attachment and optimally resistant against intelligent attack). 



5.4- Connection with results of previous studies - fraction of removed nodes 

Previous studies of network resilience, based on analysis of static topological properties 
of networks under attacks, had shown that power law networks (such as the ones 
produced by a preferential attachment mechanism) are more resistant than Poissonian 
ones against random removal of a fixed fraction of nodes (see e.g. [HI US] ) , but are very 
vulnerable against hub removal [I1I21]- In the light of our new results, one may wonder 
how the fraction / of nodes removed varies among different degree distributions, when 
considering attacks constrained by degree dependent node removal costs, with limited 
attack resources. The results of numerical explorations for C = (where optimal 
attacks will target hubs) and q — 3 are shown in Fig. [3l in the form of scatter plots 
of the relative variation AF/F — (Fbofore — raftcr)/rbcforo of the integrity measure 
under optimal hub-targeted attacks versus the fraction of sites removed (left), and 
under random attacks versus the fraction of nodes removed (right). The dotted line 
in the latter plot shows the theoretically predicted (constant) value of the relative 
variation of the process integrity measure under random attacks for C = 0. Figure |3] 
reveals that biological PINs are indeed affected by a dramatic drop in the integrity 
measure under hub removal, even for tiny fractions of removed nodes; this confirms our 
intuition that their observed resilience depends crucially on having degree statistics 
such that attack costs prevent intelligent attackers from removing significant numbers 
of hubs. The same statement is expected to apply to any random graph drawn from 
the ensemble ^, where the imposed local degrees are those of the biological networks. 
Finally, Fig. [3] also shows that, as expected, the fraction of sites removed during hub- 
targeted attacks in a Poissonian graph, with fixed attack resources and when the node 
removal cost function is monotonically increasing with fc, is considerably larger than 
in power law graphs. We conclude that the often claimed superiority of Poissonian 
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Figure 4. Left: process integrity measure r[p,g] before and after optimal 
intelligent and random attacks, when the node removal cost function is a constant, 
(/>(fc) = (k) — 1, and attack resources are constrained according to k~^ = {ij>{k))/5. 
The theoretical upper and lower bounds of 14011 are shown as dotted and dashed 
lines, respectively. Right: Scatter plot of the relative variation AF/F of the 
integrity measure versus the fraction of sites removed under intelligent attacks, 
with attack resources constrained according to k~^ = {<f>{k))/5 and with node 
attack cost function <^(fe) = (fe) — 1. Here different markers correspond to different 
network families, similar to Fig. |3] 



networks over power law graphs for hub targeted attacks is strictly a consequence 
of the decision to keep the fraction of removed sites fixed. This is consistent with 
the findings in fS], where it was argued that power law networks are no longer more 
fragile than Poissonian graphs against hub targeted attacks when one looks at the 
number of removed links, and that the efficiency of hub removal in power law graphs 
would mainly lie in the fact that this removes many more links than it would have in 
Poissonian graphs. 

In order to make contact with earlier results in literature, we consider below 
optimal attacks calculated for the constant cost function (j){k) = (f). Here the effects of 
attack costs should vanish from the problem, and our attacks should reduce to those 
where the fraction / of degrees to be removed is kept fixed. In fact, from the resource 
constraint one has / = J2kP(^)'ii^\^) ~ I/kc/). We plot in Fig. |4] (left) the integrity 
measure r[p, g] for the different networks considered so far, before and after optimal 
intelligent and random attacks, with constant cost function (/) = (k) — 1 and resources 
K~^ = (p/b. In the right panel we show a scatter plot of the relative variation of 
the process integrity measure under optimal intelligent attack versus the fraction of 
removed sites, similar to Fig. |31 The fraction of removed sites is now constant, as 
expected, and indeed equals 1/5 for our choice of the resource limit. We see that for 
the constant node removal cost function the dependence of the post-attack integrity 
measure on the degree distribution is drastically different from that in the case of 
monotonically increasing cost functions, and we retrieve the old results known from 
literature: power law networks are now more resilient than Poissonian ones against 
random attacks, but are extremely sensitive to hub removal. 
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Figure 5. Log-plot of the post-attack integrity measure for the case where 
assumed and actual node removal costs and resource constraints need not be 
identical. Our network has the size and average connectivity of the C. Elegans 
PIN. The node removal cost function is taken from a family 0f (fc) and the resource 
constraint is k" = {tj>i{k))/qi, with £ g {1,2,3,4} (see main text for further 
details). The defender assumes that I = d and chooses the associated optimally 
resistant degree distribution, whereas the actual value \s I = a, and the attacker 
bases his strategy on the latter. The defender is optimally prepared only for a = d. 



5.5. Misinformation 

We finally illustrate briefly the possible effects on the network resilience problem of 
misinformation, i.e. a situation where a network is designed to be optimally resistant 
against an optimal intelligent attack on the basis of a node removal cost function 
(j)d{k) and a resource limit k^ , but where in fact it faces an optimal intelligent attack 
constrained by an actual cost function 4>a(k) and with resource limit k~^ . Here the 
cost functions (j>e{k) and resource limits kJ are defined as follows (with integer ^ > 1): 
^i(fc) — (fc) — 1, 0£>i(fc) — fc^^^(fc— 1), and kJ — {(j)e{k))/q( with qi — 5 and qiyi — 3. 
Note that for a = 3 the intelligent attacks in fact reduce to random ones. The results 
of our numerical explorations are shown in Fig [5j For every attack, the distribution 
for which the post-attack process integrity measure is largest is indeed seen to be the 
one which is optimally resistant to the actual attack, i.e. the choice d = a (for d = 3, 4 
the optimally resistant degree distributions are very similar, and their behaviour is 
almost identical). Fig [5] suggests that, as long as the node removal cost functions are 
monotonically increasing with the node degree (i.e. for a > 2), networks which are 
optimally resistant against non-hub attacks (a = 3, 4) are reasonably resistant against 
hub-targeted attacks (a = 2), whereas networks which are well prepared against hub- 
targeted attacks behave quite poorly when subjected to non-hub attacks (a = 3, 4). In 
other words, degree statistics designed to be optimally resistant against hub removal 
appear to be quite sensitive to misinformation, whereas those optimally resistant to 
random attacks suffer less from misinformation, at least as long as the node removal 
cost function is monotonically increasing with the node degree. 
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6. Discussion 

Many research papers have been devoted recently to the resihence of networks under 
attacks. Most study resihence in terms of the behaviour of static properties of 
networks under random and intehgent removal of a fixed fraction of sites or bonds. 
Results obtained empirically [1] [25] or analytically (within mean-field and asymptotic 
approximations) [71 [2 [26l |6l [27] have shown that power law networks are more resistant 
than Poissonian ones against random attacks (see e.g. |22[|23) V but are very vulnerable 
against hub removal [I] [23]. In contrast, more recent studies [3] suggest that power 
law networks are not more fragile than Poissonian graphs against hub-targeted attacks 
when one looks at the number of removed links (as opposed to nodes) . 

In this paper we have sought to study network resilience in a more realistic setting, 
where attackers have fixed resources and where removing or disrupting a node carries a 
cost for the attacker which depends on the degree of the disrupted node. We quantify 
the resilience of the system in terms of the process for which the network acts as 
infrastructure, based on determining the critical temperatures for the onset of various 
types of global order that could be envisaged (the resulting network integrity measure 
is only weakly dependent upon the specific choices made) . This formulation also allows 
for attacks involving partial disruption of individual nodes, which would have been 
inaccessible to the techniques normally used when studying network resilience, such as 
percolation theory. We can define precisely the most damaging attack strategy, given 
knowledge of the degree sequence of a network, and for any given node removal cost 
function. In addition we could subsequently define the optimal network topology, i.e. 
the degree distribution for which the integrity of the collective process is preserved 
best when attacked by a foe who employs the most damaging attack strategy. 

A network's resilience against attacks is extremely sensitive to the dependence of 
the node removal cost function on the degree of the targeted node. This dependence 
determines the crucial outcome of the competition in such scenarios between the 
benefit and the cost of attacking high-degree nodes. If we choose a trivial constant 
cost function, we retrieve results from literature on network resilience under random 
and targeted removal of fixed fractions of sites or bonds. However, as soon as one 
chooses more realistic node removal cost functions, that increase sufficiently fast with 
the node's degree, power law networks are found to be more resistant than Poissonian 
ones, even against optimized intelligent attacks. Our results show that "modular" 
configurations with a core of nodes highly connected to each other, in a sea of (almost) 
disconnected nodes, and strongly disassortative configurations, are, depending on the 
attacker's resources, the most resistant ones against hub-targeted attacks, respectively. 
Broad distributions with fat tails are the best defence against random and low degree 
targeted attacks. We also touched briefly upon the effects of misinformation, where a 
network is designed to be optimally resistant to a certain attack, whereas it actually 
faces a different one. Results suggest that for monotonically increasing cost functions, 
degree distributions with fat tails are much less sensitive to misinformation effects. 

Upon comparing real protein interaction networks with random networks of the 
same size and average degree, we found that the attack resilience of the biological 
networks is superior to that of power law and Poissonian ones, even against optimized 
intelligent attacks. It may be that topological properties beyond the degree sequence 
play an important role here, and this deserves further investigation. In particular, one 
could calculate the integrity measure for processes supported by networks drawn from 
ensembles tailored to the production of graphs with built-in structure beyond that 
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imposed by the degree distribution, along the hnes of |2H1 H] • Another direction for 
future work may be to consider graph ensembles in which both the network topologies 
and the node removal cost functions involve hidden variables. 

Our paper emphasises the importance of distinguishing between different classes 
of network attacks on the basis of the node removal cost function and resource 
limitations imposed upon the attacker, and of studying and quantifying network 
resilience strictly within a given class of attacks. Previously proposed conclusions 
about the vulnerability of power law networks against intelligent attacks should be 
moderated in all cases where there is no compelling reason to assume that the cost to 
attackers of node removal is independent of the node degrees. 
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Most of the derivations in both appendices follow the lines of similar calculations in 
e.g. [28l [29], and we will hence be brief and highlight only crucial steps to indicate 
the changes generated by the introduction of the attack variables {i^,}- Following [T5] 
we use the property that with (k) — N~^ J^t ^i the ensemble ([4]) is identical to 



Prob(c) 



^k,k(c) 






n[^ 



c i<i 



(A.l) 
(A.2) 



Appendix A.l. Derivation of saddle-point equations 



We write the Kronecker 5s of the degree constraints in integral form, and wc introduce 
the short-hands cTi ~ {aj, . . . , af) S {—1, 1}" so that 



fA = lim lim \ log Z — log >^ / TT — -i 

^^ N^oon^o/3nN\ ^ » Z^ / llL27r 



iuiki 



lim lim-i-(logZ-log y r 17 r^e'"''=-" 



i<j 



N 



X exp 



cri...crN ' 



2N 



J2[ /dJ p(J)e/5^«'«^'^'-'^^-'('^'+"^)- 1] + OiN°)\ \ (A.3) 



We proceed by introducing for cr E {^li 1}" and f £ S the functions 



-D(^,cr|{crj,w,,e,;}) = N ^J2i^i,ii^cr.cr,e 



They are introduced via the 



substitution of integrals over appropriate ^-distributions, written in integral form: 



1 = 



d£)(^, cr)d£)(^, cr) iND{$,cr)\D(i,cr)-D(i,cr\^(Ti 



.})] 



2Tr/N 
Upon using the short hand {dDdD} = Y[^ cr ^(^' ''■)d-D(e, cr) we then obtain 



(A.4) 



Ia 



lim lim -^flogZ-log f{dDdD}e^''^^'r^(^''^)^(^''^^--^^('^)+^('°sN) 
t-,oon^o j3nN I ^ ^J^ ' 

< exp [l{k)NY, E D{e,rT)D{e',a') /dJ P( J)e^'^"''^' 



xexp[7VEp(C,fc)logE I ^e-'=-^«''^)°"'^ 



^ duj 



} 



(A.5) 

We next define z — lirajq^oo N^^ log-Z (anticipating this limit to exist), which allows 
us to evaluate / by steepest descent: 



Ja 

fn,A[ 



= lim -cxtrr^j^, /„,A [{£',£'}] 



(A.6) 



io- 
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ii' CTCT' 

+ ^p(e,fc)logE r|^ei^fe-ii)(c,<T)c-| (A.7) 

Extremization (jA.7p with respect to {_D, £)} gives the saddle-point equations: 

m, a) = m ^ Y. ^(^' '^') /d-^ P( J)e^-'«''^'^' (A.8) 

C (T' •' 

D(^,(t] = Ypitk) -^"^ ^ ^ (A.9) 

The second of these equations is simphfied using the identity 

" d^ ^i.£-ii)«,^)c-- ^ I 2n[~iD{(, ct)Y/£1 if £ > ^^^^^^ 

So, if we also re-define D{^, cr) = i(fc)F(^, cr), we arrive at 

^(^' '^) = E E ^(^' '^') /d-^ P( J)e^^«''^'^' (A.ll) 

e CT' -^ 

We note that J^e J2cr ^(Ci ^)P{^: cr) = 1 at the saddle-point. The term z ~ 
limjv_s.oo -^"^ logZ measures the number of graphs in the ensemble. It follows from 
lim^_^o(/3/) = —log 2, giving z = (fc) log(fc) — (fc) — ^j,p(fc) logfe!, and hence 

7a = - lim-l^p(e,fc)log[^F'=(e,<T)j (A.13) 

^ Ck cr 

Appendix A. 2. Replica symmetric theory 

To take the required limit n — > in our formulae we make the replica-symmetric 
(RS) ansatz. The order parameter D(^, cr) must now be invariant under all replica 
permutations, and thus have the following form: 

m.^)=j^hDi^M,^^^^^^. (A.14) 

Via equations (|A.11IA.12|) one then finds a similar structure for F{k, cr), 

F{^, (t) = fdh F{^, h) e^''I^« "° (A.15) 

and in the limit n — > 0, after some standard manipulations, a closed set of transparant 
equations for the RS order parameters D{^, h) and F{^, h): 

F{^, /i) = E /d/i'd J D{C, h')P{J)5 /i- ^atanh[tanh(/3 JCC') tanh(/3/i')] (A.16) 
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We note upon integrating and combining these equations that Jdh -F(^, h) = 
Y^^, Jdh D{C, h) = 1. This enables us to write F(C, h) = F{h\£,) with /d/i F[h\^) = 1, 
which gives immediate probabihstic interpretations of the functions F{h\S,). Upon 
eUminating I?(^, h) the RS saddle-point equations then take the new form 

Fm) = Ep(^'^)4 Ujnj) jxi[Ah,F{hM')] 

k(,' \ I J J f>^f. 

X jk- -atanh[tanh(/3J^^')tanh(/3E/i£)] (A.18) 

^ e<k 

Clearly F{h\0) = 6{h). To identify the relevant observables and calculate for 
(T e {-1,1}" the quantity P{^,k,(T) = limN^oo N'''^J2i{h,i,^k,ki^cr,cri) one uses 
the alternative form of the replica identity, viz. 



(5(^)) 



E^5(T)e-^^('^) 



= lim y 5(^1)6-'^ S^=i^('^°) (A.19) 



Upon also making the RS ansatz this results in 

Pks(C, k, cr) = pi^, k) j dh Wm, k),^-^^^^^ (A.20) 

wm,k) = j\{[dh,F{hAm[h-Y.^i] (A.21) 

•' i<k e<k 

The measure iy(/i|^, k) is the effective field distribution for those sites where {£,i, ki) = 
{^,k). We note that W{h\0,k) = S{h). With W{h) = EikPitk)W{h\^,k) we 
can write the conventional scalar order parameters m = liniTv-i-oo N~^ J2i i'^i) a-nd 
q = limjv-s-oo N^^ ^,- (cr;)^ in their familiar forms 

d/i Vl^(/i)tanh(/3/i), q ^ j dh W [h) iaW {I3h) (A.22) 

The subset of sites with {^i,ki) = (^, fc) can be regarded as sublattices in the 
sense of [30], and we can define sublattice magnetizations m(^, fc) via m{^,k) = 
Jdh W{h\£^,k)ta,nh.{Ph), such that m = X^ffe^C^' ^)"^(^' ^)- ^^ ^^^ limit T — > cx) 
(i.e. /3 — >■ 0) the only solution of (jA.18[) is as always the trivial paramagnetic (P) one: 
F(h\£^) — S{h). This is a saddle-point at any temperature, but can become unstable 
in favour of ferromagnetic (F) or spin-glass (SG) states as T is lowered. 

Appendix A. 3. Continuous phase transitions away from the paramagnetic state 

Continuous bifurcations away from the trivial state are found as usual by expanding 
(|A.18|) in moments of F{h\^), assuming the existence of a small parameter e with 
0< |e|<Cl such that Jdh h^F{h\^) = 0(6^). With some foresight we define a function 
7(^) and two |S| x |S| matrices M(^)(/3) with entries Mf^,{P), for I € {1, 2}: 

7(0 ={k)-'Y.p{tk)kik~l) (A.23) 

k 

M§{13) = 7(0) JdJ P{J) tanh^(/3JC0) (A.24) 

We define Amax(/5) as the largest eigenvalue of Af (^^(/?). If the first order to bifurcate 
away from F{h\£,) — S{h) is e^, the bifurcation is towards a state where m f^ 0, 
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i.e. describing a P-^-F transition. Upon multiplying both sides of (jA.lSp by h and 
integrating over h, the bifurcation condition for this is found to be 

P ^ F : AWJ/3) = 1 (A.25) 

If instead the first order to bifurcate is e^, the bifurcating new state has m — and 
q > 0, describing a P— >-SG transition. Upon multiplying both sides of ()A.18|) by ft,^ 
and integrating over h, the bifurcation condition for this is found to be 

P ^ SG : A(,?L(/3) = 1 (A.26) 

We focus on a specific simple bond distribution, the binary P{J) = ^(1+77)(5(J— Jq) + 
i(l — 77)(5(J+ Jo) (with Jo > 0), where the matrices M^^'>{(3) take the simple form: 

A^^V) - ?7tanh(/3JoCC')7(C') Mg)(/3) = tanh2(/3JoCC')7(C') (A.27) 

Appendix B. Equilibrium analysis for model B 

Appendix B.l. Derivation of saddle-point equations 

The calculation for coupled oscillators is initially very similar to the previous one, with 
summations replaced by integrations. The main differences start at the introduction 
of the replica-symmetry ansatz; from then onwards we have to implement appropriate 
adaptations of the calculation for XY spins in [2^ (an alternative route would be 
to adapt the cavity-based analysis in [31]). As before we write degree constraints in 
integral form, and we introduce the short-hands 0,; = [O], . . . , 0") € [— vr, tt]" so that 



X exp 



■J2^fdJ P(J)e^^«'«^I^=^°^('''°~''^°)-'("'+"^)-l] +0(A^")]} (B.l) 



.2A^' 

We next introduce for 9 E [— tt, tt]" and ^ G {0, 1} the functions D{^, 6\{6i,uJi, Ci}) = 
N~^ J2i ^^.^i^Wi 0j:]e^'"', via the substitution of functional integrals over appropriate 
(5-distributions, written in integral form. With the short hand {dDdD} = 
Y[^ D{e, 6)dD{e, 0) we then obtain an expression in the form of path integral: 

7^= lim lim-i-|logZ-log/{di?d^}e'^S.^^^(^'^)^(^'^)-^^<''^>+°('°s^) 
X exp \^{k)NY, fdede'Die,e)D{e',e') f dJ p(J)e'^^"'E„ ™«(»°-»=) 

ee' 

xe.p[NY^p{C,k)logJ\ef^^e^-'--^^i^^(^^^-]} (B.2) 

With z = limAr_>oo N~-^ logZ^v = (fc) log(fc) — (fc) — J2kPi^) logfc! (which has already 
been calculated earlier), we evaluate / by steepest descent: 

7b = lim-extr,^^,/„,B[{i^,i^}] (B.3) 

n— !-0 n x^^^s 

/n,B[...]= --^[^Y, j do b{(,,e)D{^,e) -\{k) - z 
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+ ^WE fddde'D{^,e)D{e,e') fdJ P{J)e^'^^'^^ 

+Ep(e,fc)iog r de r^e'""^'^«'^)°"'^| 
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cos(eQ-e^) 

(B.4) 



Functional variation of (IB.3|) with respect to {D, ZJ}, followed by application of (jA.lOp 
and transformation via I?(^, 9) — i(fc)F(^, 0), gives the saddle-point equations 



F{t6) = E /"de' D{^,e') /"dJP(J)e^''««'S„™^(' 



fc F''-^{^,e) 



Again ^^ Jd0 -D('C, 6)F{£^, 0) = 1 at the saddle-point, and we obtain 



fs= - hm ^ 5] P(e, fc) log d9F\^,e 



?fc 



(B.5) 
(B.6) 

(B.7) 



Appendix B.2. Replica symmetric theory 

For real-valued variables the replica-symmetric ansatz is less straightforward. 
Permutation invariance with respect to components now implies that -D(^, 6) and 
F{^, 0) are functional integrals over the space of normalized functions P : [— tt, tt] ^- IR 
(i.e. /^^d6' P{e) = 1), with functional measures VFr,[C,{P}] and Wf[£,AP)]- 



Di^,e)= f{dP}WD[t{P}]llPiOc.) 

-' a 

Fit 6) - J{dP} WF[t{P}]l[Pi9„) 



(B.8) 
(B.9) 



(we may use the same symbol P as employed to define the bond probabilities via 
P{J), the arguments will always prevent ambiguity). Insertion of (jB.8IB.9[) into the 
two equations (jB.5IB.6p then gives, in the limit ti —> and after some manipulations, 
the following closed equations for the RS measures Wd[Ci {P}] ^-^d WfK, {P}]- 



WF[t{P}] = J2j{dP'}WD[eAP'}]JdJPiJ) 



27rIoil3JCe) 



(B.IO) 



WD[t{P}]^Yl 



Pitk)k/{k) 



fc>0 [j{dP'}WF[t{P'}] 



ll[{dPe}WF[t{Pe}] 



l[s[p{e) 



Ue<kPdd) 



Ide'Ue<kPd&') 



(B.ll) 



Functional integration of both equations over P shows that /{dP}Vl/F[^, {-P}] ~ 
J2^' HdPjWoWAP}] = 1- This allows us to write WfI^IP}] = Wf[{P}\^] with 
J{dP}Wi?[{P}|C] = Ij which allows also here for probabilistic interpretations of 
the order parameters M^f[Cj {-P}]; which are now functionals acting on the space of 
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probability distributions over the interval [— 7r,7r]. Upon eliminating W£)[^, {P}] the 
RS saddle-point equations then take the following form (where 6,6' £ [— tt, tt]): 

WrRPm] = Ep(^''^)7ry /d^ PiJ) /n [{dPi}WF[{Pi}\a 

T-r r , ^ Jd9' ePJii''°^('''")lli 
xYls\P{6)- ■' ^^^ 



\.e<k 



Pim 



27rIoi/3J^^')Jde'Ue<kP^ 



(B.12) 



We observe that Tyi.[{P}|0] = Ue^iPi^) ~ (27^)"^]. To identify the physical 
meaning of our observables we define and calculate the quantity P{^, k, 6) = 
limpf^oo N~^J2i {^i,ii^k,kiS[0,6i]). Within the RS ansatz it is found to be 



PRs{i,k,e) =p(e,fc) f{dP} W[{P}\C,k]l[P{6^) 

w[{pm,k] = [n [{dPe}wA{Pim]\ U^[pi()'^-j^ 
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(B.14) 



The functional measure T4^[{P}|<^, fc] generalizes the concept of an effective field to 
an 'effective' angle distribution of those oscillators with {^i,ki) = {^,k). Note that 
W[{P}\0,k] = UeS[P{0) - (2^)-!]. With W[{P}] = j:ikPi^^k)W[{P}\C,k] we can 
write the conventional types of scalar order parameters in a compact form: 



I™ ^E (/(^')) = /{dP}w^[{p}] Tde Pid)f{6) 



N 



hm Y.{f{6,)){g{6,)) 
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{dP}w[{p}] / de* P{6)f{e) / de P{e)g{e) (b.i6) 



For T -> oo (i.e. /3 — ?> 0) the only solution of our equations is the trivial P state of fully 
random phases 6,: Wf[{P}\^] = W[{P}\£„k] = W[{P}] = Y{g6[P {6) -{2tt)^\ which 
solves out equations at any temperature, but will destabilize at some T in favour of 
ordered states with (partially) frozen relations between the phases of the oscillators. 

Appendix B. 3. Continuous phase transitions away from the incoherent state 

To find continuous bifurcations away from the incoherent (P) state one has to carry 
out a Guzai expansion '29' of the functional order parameter equations (JB.12P around 
the solution VKi?[{P}|C] ^]\g5[P{6)-{2TT)-^]. This wih involve the modified Bessel 
functions /„(z) [k]. One writes P{6) = (27r)-i+A(6') and Wpl JPW ] -> iy[{A}|5], 
with T:4^[{A}|^] = as soon as ^'^^d6 A{6) ^ and one expands (|RT^ in A(6'): 

W[{Am = ^p(e',fc)A /'dJP(J) [H \{dA,}W[{Aem'] 
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We next evaluate functional moments of both sides of this equation. If the first 
bifurcation away from the P state is of order A, we multiply by A(6) and integrate 
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(functionally) over all A, leading to an eigenvalue problem for the functions ^{(0) — 
J{dA}W[{A}\^]A{e) subject to the constraint J^^dO ^^{6) = 0: 

with 7(^) as defined in ((X23)) . The solutions are of the form 'i'^(6') = V(C)e™^, with 
m £ {1, 2,3,.. .} and with ip{^) to be solved from the eigenvalue equation 

0{A) bifurcations : ^(^ - J] { J dJ P(J)^i|^)7(^')V^(C') (B.19) 

For m — 1 the bifurcating state (F) is one where the oscillators synchronize (partly) 
to a preferred overall phase, whereas for m > 1 the transition is towards a state with 
non- uniform phase statistics but without global synchronization 29J. 

If the first bifurcation away from the P state is of order A^ rather than A, 
so /{dA}W^[{A}|^]A(6l) = 0, we muhiply (|B.17P by A(6li)A(6'2) and integrate over 
all functions A, leading to an eigenvalue problem for the function '^^{9i,02) = 
J{dA}W[{A}\i]A{e,)Aie2) subject to J^^dO, ^^ie,,e2) = /:^d02 *5(0i,02) - 0: 

*' " (B.20) 

The solutions are of the form *^(6'i,6'2) = ^(^)ei('"i»i+™282) ^ith mi^2 S 
{1, 2,3,.. .}, representing new states with 'frozen' local phase ordering but no global 
synchronization, i.e. spin-glass type states (SG), each bifurcating when 

(B.21) 

The right-hand sides of both (|B.19p and (IB.2ip vanish at /? = 0, so the transitions 
correspond to the smallest /3 such that solutions of (jB.191) and (IB.21[) exist. Hence 
we need the maxima of the right-hand sides over m and (7711,777.2)) respectively. The 
properties of the modified Bessel functions (see e.g. [IHI) ensure that these maxima 
are found for 777 = 1 and (7771,7712) = (1,1)- Finally, if we again choose the bond 
distribution P{J) = ^{l + ri)S{J~ Jq) + ^{l — ri)S{J + Jo), the bifurcation conditions 
can once more be written in the form (IA.25IA.26p . but where in the case of coupled 
oscillators the largest eigenvalues Ajnax and Amix refer to the following matrices 

Comparison with (JA.27P shows that, inasmuch as the location of the transition lines 
away from the P state is concerned, the differences between having interacting Ising 
spins or coupled oscillators on the nodes of the network are accounted for by the simple 
substitution tanh(z) — > Ii{z) / Iq{z) in the relevant remaining eigenvalue problem. 
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